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Abstract 

We derive an analytical expression for the two-gluon production in the pA (light- 
heavy) collisions, and focus specifically on the rapidity dependent part. We approx- 
imate the gauge field from the heavy target as the Color Glass Condensate which 
interacts with the light projectile whose source density allows for a perturbative 
expansion. We discuss the longitudinal correlations of produced particles. Our cal- 
culation goes in part beyond the eikonal limit for the emitted gluons so that we can 
retain the exponential terms with respect to the rapidity difference. Our expression 
can thus describe the short-range correlations as well as the long-range ones for 
which our formula is reduced to the known expression. In a special case of two high- 
Pt gluons in the back-to-back kinematics we find that dependence on the rapidity 
separation is only moderate even in the diagrammatically connected part. 



1 Introduction 



The gluon distribution function in hadrons grows up at high collision energy 
and thus small Bjorken's x, which makes it indispensable to elaborate a re- 
summation over the multiple scattering to abundant soft gluons. When the 
transverse distribution of partons is saturated eventually at small x and small 
Q 2 [1], the classical approximation is to be a good description of the hadron 
wave-function dominated by coherent small-x partons. In fact, a shift in the 
perturbative QCD (Quantum Chromodynamics) vacuum by such a classical 
field is an efficient prescription to take account of dense gluons within the 
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eikonal approximation. This resummation scheme has been well founded in 
the non-linear but still perturbative regime embodied by the notion of the 
Color Glass Condensate (CGC) [2,3]. 

In this paper we are interested in the inclusive one- and two-gluon produc- 
tion cross section in the CGC formalism. It is not quite viable to solve the 
non-linear problem fully in the case that the target and the projectile both 
have dense gluon contents. In contrast, the analytical calculation is feasi- 
ble if we deal with the asymmetric case of a light projectile scattering off 
a heavy target [4], as is the case in the pA collision or in the forward and 
backward rapidity regions in the AA collision. We remark that the inclu- 
sive multi-particle production has been discussed in several different frame- 
works [5,6,7,8,9,10,11,12,13,14,15,16]. The present work aims to elucidate the 
production rate of two gluons having two distinct rapidities y\ and y 2 in the 
asymmetric collisions. We have to limit the separation \y x — y 2 \ parametri- 
cally less than l/a B since we do not take account of quantum evolution which 
is incorporated by the radiative corrections between two separate gluons. It 
is maybe worth noting that, nevertheless, we can describe another aspect of 
quantum evolution through the hadron wave-function, that is, the evolution 
of the saturation scale [17], though it is non-trivial whether and how these 
quantum corrections are to be factorized [16]. 

We shall accomplish the whole calculation of the two-gluon production along 
the line of Refs. [9,18,19] and thus our results should be an extension of that 
given in Ref. [8] where yi ^> y 2 is assumed for simplicity. Our results are in fact 
consistent with Ref. [10] in the eikonal limit that emitted gluons are energetic 
enough to overwhelm the transverse recoil, while they possess intricate extra 
contributions depending on Ay = y± — y 2 with the eikonal limit relaxed. We 
will later discuss what the missing part was that causes the discrepancy. 

It should be an important problem to clarify the particle correlation reflecting 
the parton saturation effect. This is because the CGC configuration naturally 
leads to strong chromo-electric and chromo-magnetic fields which extend be- 
tween the target and the projectile [20,21]. It would be conceivable that the 
longitudinal correlation provides us with information on the initial dense state 
to be regarded as CGC matter. 

Moreover, in Ref. [22], the jet azimuthal back-to-back correlation has been 
considered, which is also discussed in Refs. [10,14]. Also, from different inter- 
est, the recent experimental data at RHIC (Relativistic Heavy-Ion Collider) 
exhibits the wide rapidity correlation in the near side of the trigger and as- 
sociated jets [23,24,25,26], which could have an influential remnant from the 
initial state and the strong color field in it [27,28,29]. In a slightly different 
context, in Refs. [30,31], the forward-backward correlation has been discussed 
which seems to have a direct link to the longitudinal field presumably existing 
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in CGC matter. Such a long-range rapidity correlations have been found in 
the recent experiment at RHIC [32,33]. All these preceding efforts are decent 
to motivate us to carry out the analytical QCD evaluation in the same pro- 
cedure as Ref. [9] where only minimal assumptions under theoretical control 
hide behind modeling. 

In the subsequent sections we shall elucidate the two-gluon correlation starting 
with the established computation for the inclusive one-gluon cross section for 
self-contained discussions. 



2 Gluon Production 

We here make a quick review on the formulation to evaluate the gluon produc- 
tion amplitude, which should be useful for clarity of our convention adopted 
in the present paper. 

2. 1 One-gluon production 

The production amplitude of one on-shell gluon with three-momentum p, he- 
licity A, and color a is expressed as the amplitude from the initial state \Qm) 
to the one-gluon asymptotic state \p, A, a). Using the Lehmann-Symanzik- 
Zimmermann (LSZ) reduction formula [34], we can generally write this am- 
plitude in a form of 

(p;\,a\n in ) =ieW(p) J d 4 xe ip ' x D x (n ont \A^x)\Q in ) , (1) 

where A^(x) denotes the gluon field operator. We note that we should take 
the on-shell limit po —>■ E p = \p\ in the right-hand side because the emitted 
gluon is real hence on-shell. The phase space integral and the summation 
over physical degrees of freedom lead us to the following expression for the 
multiplicity of emitted gluons per rapidity 0; 




Precisely speaking this expression is not the multiplicity but the probability for 
the one-particle emission. Nevertheless, this simply coincides with the multiplicity 
in the leading order of g in the present pA problem. We do not distinguish the 
multiplicity and the one-particle emission probability here for this reason. In the 
case of the AA collision, however, all multiple-scattering processes contribute to the 
multiplicity [35]. 
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Here we have introduced the momentum rapidity variable as y = \ \n[p + /p ] 
where p ± = (p° ±p z ) / a/2 which can be expressed as e ±y / \/2 by using the 
on-shell condition, i.e. 2p + p~~ = |p_|_| 2 . It is then straightforward to confirm that 
the phase space integral is rewritten as d 3 p/[(2n) 3 2E p } = dyd 2 p ± /[2(2n) 3 ]. In 
Eq. (2) the summation with respect to A should run over only the transverse 
components corresponding to physical degrees of freedom, and (• • • ) denotes 
an average over the color configuration of sources. 

It is convenient to adopt the following convention to define the physical po- 
larization in the light-cone coordinates [36]. Let us particularly consider the 
light-cone gauge defined by 

n^ = A+ = 0, (3) 

where n M = 8 + . Then, as beautifully demonstrated in Refs. [18,37], we can take 
advantage of this gauge choice to reduce computational cost significantly. The 
light-cone gauge condition leads to e^ A ^ + = and there are not four but three 
independent polarization vectors, two of whose linear combinations are phys- 
ical. The following projection operator enables us to identify the transverse 
parts simply with A = 1,2, that is, 

e£\p) = -D» x (p) , (4) 

where 

D^(p) = -g^ + — - P n^n v . (5) 

n-p (n ■ p) 1 

The above tensor structure appears also in the gluon propagator in the light- 
cone gauge. It is important to note that these polarization vectors are doubly 
transverse to both the momentum and gauge directions, i.e. 

p^(p) = 0, n"eW(p) = 0. (6) 

Moreover, using the explicit form of Eqs. (4) and (5), we can show that the 
sum over physical polarization amounts to 

E #>{p)4 x >{p) = D lw {p). (7) 

A=l,2 

Since the above is multiplied by A^ 1 and A v after all (see Eq. (1)), only a 
quantity of A ll D llv A v is necessary to reach the final expression. This quantity 
is reduced to as simple as A 1 A 1 owing to the gauge condition (3). Therefore, 
we only have to retain a 5ij term out of D^p) in the polarization summation 
to estimate the one-gluon production. Consequently the gluon multiplicity is 
expressed as [9,18,38,39,40,41] 

(^)-5/&<W<*m*<p>>. («) 
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where we have denoted the vacuum expectation value of the gauge field, i.e. 
background gauge field, as A%(x) = (Q\A£(x)\Q) . Because p 2 in Eq. (8) is a 
four- momentum squared, only the singular contribution proportional to 1/p 2 
in A l a (p) remains non-vanishing in the on-shell limit that we should take in 
the end. We note that we define the Fourier transformation in our convention 
as 

j£{p) = fd'r^iA'Jr). (9) 
Generally this is a complex- valued function and satisfies A l *{p) = A l a (— p). 



2.2 Two-gluon production 

It is a straightforward generalization to continue on considering not only single 
but also the multiple gluon production in this way. We will focus particularly 
on the two-gluon case in this work, though it is in principle possible to ac- 
commodate more gluons. In the same way as previously explained, we can 
write down the production amplitude for two gluons with three-momenta p, 
q, helicity A, a, and color a, b, respectively. The LSZ reduction formula yields 
the amplitude as 

(p,\,a;q,a,b\Q iQ ) 

= i C W(p) J dW-a i$\q) J d 4 ye^D y (Q out |T[^(x)^(y)]|Q in ) . 

(10) 

Here T represents the time ordering as usual. From this amplitude we can 
reach an expression for the two-gluon correlation in the leading order of g for 
one emitted with rapidity y\ and the other emitted with rapidity y 2 ; 

The first term corresponds to a long-range correlation part which is not sup- 
pressed for large \y x — y 2 \. This term is important to consider the azimuthal 
angle correlation, which arises through the average over the CGC configura- 
tions [28]. The second is an interference term, and the third is a short-range 
correlation part which depends on Ay as diminishing as |Ay| goes large. To 
extract the connected contribution henceforth, we have defined the propagator 
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(A) (B) (C) 




Fig. 1. Diagrammatic representation of three distinct contributions to the two-gluon 
production: (A) the disconnected O(pp) diagram, (B) the interference gO(p^) one, 
and (C) the connected g 2 0{p^) one. 

and its Fourier transform as follows; 

Q%{x,y) = (n out \T[A^x)Al(y)}\Q in ) - AZ(x)A» b (y) , (12) 
QZ(p, q) = J d 4 x d% e^-y g%(x, y) . (13) 

Just for clarity we remark that GH\q,p) literally means g*^(p,q). The cen- 
tral discussion in what follows lies in analytical evaluation of the background 
gauge field in the CGC formalism and the associated propagator on top of the 
background. 

We can pictorially visualize these three terms enumerated in Eq. (11) in a way 
as sketched in Fig. 1. The first (disconnected), the second (interference), and 
the third (connected) terms of Eq. (11) correspond to Figs. 1 (A), (B), and 
(C), respectively. The cross represents the interaction with the light projec- 
tile whose intrinsic color distribution is regarded as small perturbation, i.e. 
p p = O(g). As a result of this setting (A), (B), and (C) have the same or- 
der; (A)= G(p 4 ), (B)= g O(fi), (C)= g 2 0(pl), and they are all O(g'). We 
can drop (B) under an assumption of the Gaussian average over p p . Although 
the long-range correlation stemming from (A), which is predominant in the 
AA collision, is physically interesting [28,30,35], we will not consider (A) and 
concentrate in (C) only that is the contribution carrying the short-range cor- 
relations like e ±Ay . In fact, what we will pursue in this work is to retrieve the 
neglected rapidity dependence in the two-gluon production with emphasis on 
the analytical formulation. We will report the full correlation including (A) 
and phenomenological implication in a separate publication [42]. 

We shall comment that the leading-order piece with the rapidity dependence 
is identified differently according to whether it is the pA or the AA collision. 
The rapidity dependence in the pA case comes from the connected diagram 
in the leading order, while the AA collision has no rapidity dependence be- 
cause the leading-order diagram is not the connected but the boost-invariant 
disconnected one. Hence, the rapidity dependence in the AA collision starts 
from the next-leading order where the logarithmic terms in rapidity appear. 
It is possible to resum all the leading-logarithmic terms only to renormalize 
the weight function of the color source distribution, resulting in the evolution 
by means of the JIMWLK equation [16]. The same treatment for the leading- 
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logarithmic terms could be applicable to the disconnected diagram in the pA 
case. 



3 Background Gauge Field with Dense-Dilute Sources 

We shall adopt the CGC picture to proceed to concrete calculations. The back- 
ground field originates from the dense color distribution inside a fast-moving 
hadron. In a classical approximation in which quantum corrections (i.e. the 
wave-function renormalization to the field expectation value) are neglected, 
A^{x) is simply a solution to the classical Yang-Mills equations of motion 
with color source pa(x±). This type of the one-source problem is solvable and 
the explicit form of A fl {x) as a function of Pa(x±) is known [43]. 

In this paper, as we mentioned in Introduction, we shall consider a situation 
with two color sources corresponding to asymmetric collisions; one source dis- 
tribution is dense from a heavy nuclei which gives the color current J~ = 
5(x + )pa(x±) and the other is dilute from a proton for instance, which gives 
J + = 5(x )p p (xi_). We treat p p (x±) as small perturbation, where we should 
note that the classical approximation still works because of the presence of 
large p A (x±). 

Here, these naive color currents do not satisfy the covariant conservation. We 
thus need to augment these definitions with appropriate gauge rotations as 

Ja = 8(x + )V ab (x)p Ab (x ± ) , J+ = S{x~)Uab(x)pp b (x±) , (14) 

where V(x) = V x - exp[ig J x dz~A + ] and U(x) = V x + exp[ig J x+ dz + A~] with 
V x ± being the path-ordering with respect to x ± . The color matrices A^ are 
the full gauge fields (including both classical and quantum parts) in the color 
adjoint representation. The covariant conservation is then manifest because 
D+V(x) = and D~U(x) = 0. 

The analytical calculation makes sense when either Pa{ x ±) or Pp{ x ±) is small 
compared to the gluon transverse momentum (i.e. saturation scale). We can 
then perform the perturbative expansion in terms of small p p (xj_). The gauge 
field associated with such an asymmetric collision can have an expansion like 

where .4^ and are the zeroth and first order terms with respect to p p (x±), 
respectively. One can in principle continue the expansion up to arbitrary order, 
but what we need to retain in this work is only up to the first order. In this 
formalism of the pA problem the gauge fields become time independent after 
interaction at the collision point, as we will see shortly. The LSZ formula 
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means that the residue of the massless pole corresponds to the probability of 
the gluon associated with the gauge fields in the asymptotic state. This makes 
a contrast to the approximate prescription often used in the AA case; the 
field amplitude for each momentum mode translates into the particle number, 
which could depend on time [44]. It is not easy how the LSZ formula works 
for the AA case [35], in particular in the expanding geometry, while there is 
no ambiguity in the pA collision. 



3.1 Gauge choice 



It is practically crucial to take an appropriate gauge choice to make the whole 
computational procedure transparent. In what follows, we will make use of the 
light-cone (LC) gauge A + = with the background field Afo which is usually 
referred to as a solution in the covariant (COV) gauge, 



A, 



(o) 



0, A, 



(o) 



A\ 0) = . 



(16) 



It should be mentioned that V(x) = 1 in this A + = gauge so that no compli- 
cation enters involving p\(x±), which is a great advantage in this technique. 
In contrast to that, even in the leading order in p p (x±) at the classical level, 
we have to keep U(x) which stems from nonzero A^ y Because A^ has 5(x + ) 
in it, the x + dependence of U(x) is given as 



U(x) 



1 for x + < , 
U(x±) = V x + exp 



r0+ 

ig dx + A( 0] 

J 



for x + > + . 



(17) 



We will use the same notation for U(x) and U(x±) hoping that no confu- 
sion would arise from this. We remark that the integration with respect to 
x + in the exponential needs a careful treatment with appropriate regulariza- 
tion [9,45,46]. 

To reiterate explicitly, our choice (that was first adopted in Ref. [18]) is 



LC for p 



A" 



COV for A 



(18) 



The essential point in this trick invented in Ref. [18] is that the classical 
solution in the COV gauge for left-moving Pa( x ±) is consistent with the LC 
gauge for right-moving p p (x±). Actually, since A^ = 0, the solution (16) 
obtained in the COV gauge for p\(x±) may well be regarded as a solution 
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in the LC gauge for p p (x±). [Usually A^ = is assumed in addition to the 
gauge fixing condition.] 

The benefits from the choice (18) are twofold: One is understood from Eq. (7), 
that is, the light-cone gauge significantly simplifies the summation over phys- 
ical degrees of freedom. The other is the fact that only one component of the 
background field (16) is non- vanishing and the remaining A7 \ takes a finite 
value at x + = only. Therefore the gluon propagation receives no effect from 
the background field except at x + = where the eikonal phase and color ro- 
tation are induced. This nice feature tremendously reduces the computational 
labor. 



3.2 Classical solution 



We are evaluating A^ that is necessary not only for the one-gluon production 
but also for the two-gluon production. From the equations of motion, D p F liV = 
J", in the light-cone gauge we have a set of equations to determine Afo, that 
is, 



- d +2 A m + d+d'Alv = J + , 



d+D m A m - d\A m + D^A^ - 2ig[d i A« ) ,A 



(0) 1 



2d + D^A l {1) - d + d l A m - d\A l {1) + &&A\ 



0, 



0, 



(19) 
(20) 
(21) 



for v = respectively. It is notable that J disappears in the second 

equation. Using Eq. (19) we can simplify Eq. (20) into 



OA (1) -2i4-4 (0) ,d+-4 (1) l -2ig\d i A (0) ,A\ 1 



l (o) 



(i) 



l (0)'^(l) 



0. 



(22) 



Here we used D^J + = 0. In the same way Eq. (21) becomes as simple as 



UA\ l) -2ig A^d+A}, 



l (0)' 



(23) 



It is not difficult to solve Eq. (23) for ^4^ with the initial boundary condition, 



A\ 1) (x + <0,x ,x±) = 6(x )-jp-Pp(x ± ) , 



(24) 



which is nothing but a retarded solution in the presence of one right-moving 
source before the collision. 



Let us then examine the connection condition across the singularity located 
at x + = 0. By integrating Eq. (23) with respect to x + from - to + , and 
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picking up a singularity contained in .4.^ (see Eq. (16)), we can acquire the 
boundary condition, 

jf^dx + £> ( o ) ^ 1) = 0. (25) 

Here we have dropped the overall d + which is irrelevant. We immediately see 
the solution to the above equation to be 



A: 



(i) 



x 



= + ,x , x±) = U(x±)A\^(x + = ,x ,x±), 



(26) 



remembering that D~U(x) = 0. In the forward light-cone region, i.e. x + > 0, 
there is no background field and thus the gluon propagation is expressed by 
the retarded free propagator A R defined by □ x A R (x) = i5^(x) leading to 



A R (x) 



—l 



d 4 p 

(2tt) 4 p 2 + ip+e 



(27) 



Here e — > + is understood as usual. It is easy to make sure that A R is nonzero 
only for x + > 0. It is straightforward to confirm from the above definition a 
useful relation, 



d + A R (x 



=0+ 



X 



(28) 



Then the solution to Eq. (23) at nonzero x + is generally to be expressed as 



- -4(i<) + -4(i>) 



(29) 



with 



A (1<) (x) = 9(-x+)9(x-)^ rPp (x ± ), 



A {1>) (x) 



y+>0 



d 4 y -i A R (x - y) 



Qi 



+ / dy~d 2 y ± \-2id+A R (x -y)]A {1) (y) 



=0+ 



(30) 



(31) 



Here «4(i<) is the past solution coming from the source for x + < 0. The first 
term in «4/ 1> - ) has the same origin but it involves a color rotation by U(x±) 
additionally for x + > 0. In the presence of the second term «4( 1> ) satisfies 
the correct boundary condition (26) in the x + — > + limit, which is obvious 
from Eq. (28). Also we can readily find that this second term produces no 
contribution for x + > when acted by the d'Alembertian □. From these 
expressions we can operate the d'Alembertian to find 



dA {1<) = -25(x+)5(x-)p Pp (x ± ) - e{-x + )e{x-)d'p p { XA _) , 



(32) 



d 



nA (1>) = 25(x + )5(x-)U(x ± )-^ Pp (x ± ) - e{x + )6(x-)d l U(x ± )p p (x ± ) 



d 



(33) 
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We can take the Fourier transform and reach 



-P 2 ^(1<)(P) = -¥ 



_(p+ +ie)(p- -ie) 



P P (P±) 



pi 



d 2 fc 



(2tt) s 



(p+ + ie)(p- + ie) 



-2jfej 



C/(fc 2± ) 



1-L, 



K l± 



(34) 
(35) 



where fcu_ + A;2± = p±. We can compute the gluon multiplicity by substituting 
the sum of the above expressions into Eq. (8) and taking an ensemble average 
over the p p (xj_) and Pa(x±) distributions. We will not do that because the 
final results are known [9,39] and the physical quantity of our current interest 
is the two-gluon correlation. 

Although it is not relevant to the one-gluon case, the two-gluon production 
requires the functional form of A7^ for the present purpose. This is actually 
where our results come to differ from Ref . [10] . There is no source contribution 
oc J M in A^y From Eq. (22) the boundary condition is inferred from 



r0+ 

L 



dx + [ D [0) A {1) - ig 



o, 



(36) 



where we have dropped the overall d + in the same way as previously. The 
above leads to the connection condition as 



A {l) {x 



+ ,x ,x±) = U(x±)A 



(i) 



,x ,x±) 



+ 



&U(xi_)] q+A\ 1) (x + = 0-, x~, x ± ) . 



(37) 



The first term in the right-hand side is zero in fact because A7 1 y{x + = 0~) = 0. 
Nevertheless, we write it to avoid confusion in order to define a convenient 
notation in later discussions. The solution for arbitrary x + is, therefore, 
given as 



A (1< y{x) 







A a>) (x) = J dy-d 2 y ± [-2id+A R (x - y)]A( iy (y) 



y+=0+ 



(38) 
(39) 



The representation in momentum space is thus A^yip) = and 



-p 2 A (1> y[p) = -i 



' d 2 fc i± -2k 1± ■ k 2 ± 



(2k) 



p + + ie 



U{k 2± ) 



P P (fci±) 
k 2 



(40) 



where ki± + k 2 ± = p± again. We note that the right-hand side has no depen- 
dence on p~ because QA^x is proportional to 5(x + ). For consistency check 
one can readily see that »4(i>) and -4(i>) obtained above certainly satisfy the 
rest of the equations of motion, namely, Eq. (19) as they should. 
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Later on, for notational simplicity, we will introduce and C M to express 
the connection condition and the solution, respectively. 



4 Propagator with the Background Gauge Field 



It is the background propagator that is necessary for the two-gluon correlation. 
We will calculate it here by means of the expansion in terms of small A^ oc p p . 
The full background propagator without expansion is defined by 

D x D x g, u - D,D U - 2igF^]g™(x,y) = g^(x - y) , (41) 

where both A7 \ and A^ enter the background field. We will use the adjoint 
matrix notation from now on, namely, a color matrix is understood as 

(A»U = {T%KU = if acb A» , (42) 

in the color adjoint representation. Then, the quantity in the square brackets 
is decomposed into the zeroth-order and first-order parts with respect to A^y 
It follows that we have 



D x D x9flu - D^D V - 2\gF t 



= Df 0) D m g^ - D (0) ^D {0)v - 2igF {0 )^u + WIV + 0(f$) , 
where we have defined the vertex matrix of order 0{p v ) as follows; 

WIV = -i(?| [(dxAfa) + 2A x 1) d x ]g^ + 2(d„A {1> ) - 2{d v A {1) ^) } . 



(43) 



(44) 



We remark that we have dropped terms irrelevant to the gluon propagator, 
that means, terms which vanish when sandwiched by the zeroth-order prop- 
agator with transversality and also with the upper + components vanishing 
due to gauge fixing. We shall define the zeroth-order background propagator 
as 

- i[D x 0) D m g, v - D (0)tl D (0)u - 2igF i0)lu ,]gff ) (x,y) = g^\x - y) , (45) 

and then we can express the full propagator in Eq. (41) in a way as Q = 
+ ^r)~\ that yields an expansion, 

g^(x,y)=g^(x,y)- J d i zg^(x,z)5T Xa (z)g^(z,y) + 0(pl). (46) 

At this stage it is tangible how we dropped terms from Eq. (43) to Eq. (44); 
irrelevant terms in ST^ U , in fact, arise from the transverse properties, i.e., 
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(A) 


^ ► 


X / 






\Pa 


y 






Fig. 2. Some typical examples for Q^.(x,y): (A) The propagation is free from the 
source singularity in the case for (x + < 0, y + < 0) or (x + > 0, y + > 0). (B) The 
singularity at x + = affects the propagation which occurs in the case for (x + < 0, 
y + > 0) or (x + > 0, y + < 0). 

n ^(o) = an d ^^(o) = ( see a ^ so Eq. (6)). In momentum space we can write 
the vertex as 

5T^(p, q) = iM(i)(P + Q)^Mp, q) (47) 

with 

r^Go, ?) = (p - q)x g^v + ^%g v \ - 2p„g\p ■ (48) 

It should be mentioned that the above vertex is slightly different from a stan- 
dard QCD three-point vertex. It might seem that the momentum conservation 
does not hold. This is because q u out of — 2p u — q v vanishes from q u G^(q) = 
and p^ out of 2g M + p^ vanishes from Q^{j))p^ = 0. All the terms involving 
disappear due to gauge fixing. We shall then concretely compute Eq. (46) 
in the subsequent subsections below. 



4.1 G^(x,y) for (x + > 0, y + > 0) or (x+ < 0, y + < 0) 



Now let us calculate the propagator in the region without crossing the source 
singularity, that is, the region with x + > and y + > or the region with x + < 
and y + < as depicted in Fig. 2 (A). There, the background propagator is 
just the free one which is unity in color space given by 



— i 



(d x d x )g, u - (dM] A™(x,y) = g"5^(x - y) . 



(49) 



In the momentum representation, we can write the time-ordered propagator 
in a concise form as 



A^(p, g )^(27r) 4 ^)(p + g)^-^ 

p z + le 



(50) 



using D^ip) defined in Eq. (5). The standard e prescription enables us to de- 
fine the retarded (x + >y + ) and advanced (x + <y + ) propagators, i.e., A^(p, q) 
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and A^(p, q) by replacement of the denominator by p 2 + \ep + and p 2 — iep + , 
respectively. 



4.2 GZ(x,y) for (x + < < y + ) or (y + < < x + ) 



As a preparation to proceed to the case with a singularity at x + = on the 
propagating path, we shall introduce the following notation. We will express 
the boundary conditions, Eq. (26) and Eq. (37), into a single form of 



A»(x + = 0+) 



x 



o~ 



(51) 



where the gauge rotation associated with the singularity is 

M) = 5)U(x ± ) , MZ = U(x ± ) , M~ = &U(x ± )±; , 
M' M'. 0. M~ = , M+ = . 

This representation is useful in the following formulation. 



(52) 



If x + < < y + or y + < < x + then the propagator passes over the source 
singularity. We show an example in Fig. 2 (B). We will derive the appropriate 
boundary condition at the singularity here. To do so, we shall pick up only d~ 
and S(x + ) terms from the equations of motion. Then we consider the connec- 
tion condition for each polarization component \i for the case of y + < < x + 
first. From the boundary condition for the gauge field (51) it is straightforward 
to write the boundary condition down for the propagator as well in a similar 
form as 



=0+ 



{x,y) 



(53) 



We note that the free propagator appears in the right-hand side for x + = CT 
and y + < 0. Using the same technique as in the previous section where we 
calculated the gauge field, we can express the propagator for y + < < x + 
with satisfying the above boundary condition as 



z+=0 



G&(x,y) = J dz~d 2 z ± [-2id+Ag(x,z)] 9xaG^(z,y) 

= -J dz'd 2 z ± [-2^Ag(x,z)]u(z ± )A»(z,y) 

In the momentum representation the above translates into 

S&M = J d'xd' y <P* +i *-ygg ) (x,y)e(-y + ) 



(54) 



(55) 
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X 




Fig. 3. Typical diagrams representing the first-order propagator which contributes to 
the two-gluon production amplitude: (A) One gluon splits into two after interaction 
with the target. (B) Two gluons receive interaction with the target. 



Here we note that the propagator poles are located only on the lower (or 
upper) half plane in terms of the complex p~ (or q~ respectively) variable, 
which reflects the boundary condition, y + < < x + . 

Similarly, the propagator for x + < < y + is to be written as 

S&M = (2^++ ,+) %+) ^^M_ u\- Pl - , J , (56) 

in momentum space. We note that the color rotation takes place in the opposite 
orientation (i.e. U^ 1 = W) to the previous case because the gluon propagation 
penetrates into the source from the opposite side. 



4-3 Propagator and amplitude 



We now get ready to proceed to evaluating the expanded propagator (46). The 
most general form is, however, not necessary for the two-gluon production 
amplitude. In coordinate space only the late time behavior at x + ^> and 
y + ^> is relevant to the production process. We can actually confirm that the 
production amplitude does not have any contribution from the region where 
x + and y + are negative. In the following argument, thus, we shall restrict the 
time arguments as positive x + and y + only. 

The production process can be separated into two distinct types according to 
the location of the three-point vertex. Figure 3 is the diagrammatic represen- 
tation for them. We can distinguish (A) and (B) out of Eq. (46) by the sign 
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of z + so that we can decompose it into 

G^(x, y) - g$(x, y)= - f d 4 z A^(x, z) 5F {>)Xa (z) A™(z, y) 

r ( 57 ) 
- J d 4 zg^ ) (x,z)5r {<)Xa (z)G^ ) (z,y), 

where we can forget about Q^(x,y) for our purpose. In the above we have 
defined, 

5r {>)Xa (z) = 8Y Xa {z)e{z + ) , 5T (<)Xa (z) ee 8Y Xa {z)e{-z + ) . (58) 

This relatively simple form of Eq. (57) is, as we stated above, valid under the 
condition, x + > and y + > 0. The first and second lines correspond to (A) 
and (B) of Fig. 3, respectively. 



44 Evaluating A mA 5T (>)Xa A av 



First we shall explicitly evaluate the contribution at z + > shown in Fig. 3 (A). 
Since the expression in momentum space is necessary to acquire the amplitude, 
let us take the Fourier transform, 

J d 4 xd 4 yJ p - x+i ™9(x + )9(y + ) J d A zA»\x,z) 5T (>)Xa (z) A™(z,y) 

D» x (p)8r {>)Xa (p,q)D™(q) ( 59 ) 



= -6(p+) 6(q + ) 



(p 2 + ie)(q 2 + ie) 



where we have used p + > and q + > in the denominator. From Eq. (10) 
the matrix element (p, A, a; q, a, b\Q{ n ) turns out to have a contribution as 

(p,\,a;q,a,b\a m ) (>) = (p) (q) SF^p, q) . (60) 

Using Eqs. (10), (35), (40), and (47) we can figure its explicit form out as 

(p, A,a;<2,<7,&|fi in )(>) 

= 9$>{p)$\q) J ^T^(p,q;k 1± ,k 2± ) 



(61) 



where fcu_+ fc2_i_ = Pj_ + Qj_ and we have used the following function; 

Jf >) (p, g ;fc 1± ,fc 21 )= (p + g)2 + i(p++g+)e C (p + ?;fci±,fcn). (62) 
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Here, we have defined 



C + (p,k 1± ,k 2± ) = 0, 

-2k 1± ■ k 2 ± 



C (p,k 1± ,k 2 ±) 



p+ + ie ' (63) 



C\p, k 1± , k 2± ) = ^ . - 2k\ , 

(P+ + le) (p + le) 

so that we can express the gauge field after the collision in a concise form as 
-P 2 ^»(P) = "i J^m P ;k 1± ,k 2± )U(k 2 y-^. (64) 



4-5 Evaluating Q'ts ST^ X a Q\ 



In the same way we can continue the calculation by taking the Fourier trans- 
formation to obtain, 

J d 4 xd 4 ye i ^ye(x + )e(y + ) J d 4 zg£ ) (x,z)8r {<)Xa (z)g^ ) (z,y) 
d 4 fci d 4 A; 2 



/ 



(2tt) 4 (2tt) 



By substituting the explicit forms given by Eqs. (55) and (56) we can write 
the above down explicitly as 



0(p+) %+) /■ d 4 ^ d 4 A; 2 



f^^**)6(p+-kt)(2K)6(-k} + q +) 



p 2 +iep+ q 2 +ieq+J (2tt) 4 (2vr) 4 v 7 v ^ i/v y v 2 ' * y (66) 
2ktD^(p)D*(k 1 ) 2k+ D°P{k 2 )D${q) 



x 



(fci + ie) + ie) 

x C/(p ± - feu.) 5Y {<)Xa (k u k 2 ) U\q ± - k 2± ) . 

The poles in ST^xa{ki, k 2 ) in the complex fcf and A;^ plane are found on 
the upper half plane. We can then perform the integration with respect to fcf 
and k 2 avoiding the pole in 5T( < )\ (T (ki, k 2 ), so that we pick up the contribu- 
tions -{2m) (2kt)- 1 S(k^ - k 2 1± /(2kt)) from \jk\ and -(2m) (2k 2 h )~ 1 5(k 2 - 
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kl ± /(2k 2 )) from l/k 2 , which leads to 



d 4 h d 4 £: 2 
(2tt) 4 (2vr) 4 ~ 1 ; 



6(p + ) 9(q + ) f d 2 k 1± d 2 fc 2 



p 2 +iep+ q 2 +ieq+J (2tt) 2 (2tt) 2 
xU(p ± - k 1± ) 5T (<)Xa (k l , k 2 ) U^(k 2± - q ± ) , 



D^(p)D*(k 1 )D°P(k 2 )DZ(q) 



(67) 



after the integration with respect to kf and k 2 . Here the delta functions 
impose that k\ — k\ — 0, kf = p + , and k 2 = q + . From Eq. (10) the matrix 
element from this contribution becomes 



/ \ no \ (A) / \ (<r) / \ /■ d 2 fci ± d 2 fc 2 ± 

(p,\,a;q,a,b\n in ) {<) = - e ^(p)eJ, >{q) J 



D^' (ki)D clv (k 2 ) 



x U(p ± - k 1± ) Sr^ X 'a'(ki, k 2 ) U\k 2L - q ± ) , (68) 



which simplifies further to be 



(2tt) 2 (2tt) 2 



x 



(feu. + k 2± ) z 



ab (69) 
T^(p,q; k 1± ,k 2± ), 



where we have defined 



T^(p,q;k 1± ,k 2± ) 



Here we recall that k l 
condition. 



L^Orw^i, k 2 )D°'\k 2 ) (h + k 2 y 



(70) 



(p+ + <?+ + ie) (fcf/ + k 2 - ie) 
= k\ ± /2p + and A;^ = fe 2 i/2g + from the on-shell 



For later convenience we shall shift the momenta as 

fell -> fc'n = fell + fc 2 _L , k 2± -> fc'_L = Pj. - fc'n + fe 2 _L , (71) 

and omit the prime on the momenta afterward. Then, the above amplitude is 
rewritten into 



(p,\,a;q,a,b\a m ) {<) = ge^(p) e<?\q) J d ^ ^ ^ 

x ^<) (p> s; p± - fc ^> + fc ± - P±) 



U(k ± )^^-U\k ± -k 2± ) 
fe 1 _l_ 



o,6 



(72) 
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Fig. 4. Momentum assignment for (A) {p,q\£kn)(>) anci (B) {p,q\Qm)(<) respec- 
tively. 

4-6 Total amplitude 



After all, the sum of Eqs. (61) and (72) amounts to the total amplitude, 
(p, A, a; q, a, b\Q in ) =ge^ } (p) e\ > (q) J — — 



X 



T^(p,q;k l± ,k 2± )T d A U dc (k 2± ) + J 



(2tt)' 
d 2 k ± 
(2nf 



(73) 



xT ( ^(p,g;p ± -fc ± , k 1± +k ± - Pl _)U{k ± )nU\k ± -k 2± ) 



lib 



Figure 4 shows the momentum assignment overlaid on Fig. 3. In the above 
expression the SU(A r c ) algebra in the adjoint representation is denoted by Ta- 
Let us now confirm that this results are certainly reduced to zero in the limit 
of vanishing target source, that is, U(k±) — > (2ir) 2 8( 2 \k±), for consistency 
check. 



Then, the quantity inside the square brackets is 

T A [ T (>) (p,q\P±+Q±,0±)+ (p, q;p ± ,q ± ) 
After some calculations we can find that 



(74) 



T (>) (P, P±+q±, 0j.) = -T^ (p, q; p ± , q ± ) 

D^(p)T Xi(7 {p,q)D^(q)(f+q i ) 



(75) 



(p + + q + )(p-+ q ) 
which makes Eq. (74) canceled as anticipated. 

For later convenience let us write Eq. (73) in a slightly different way. Noting 
that, 

|^ U(k ± )TZtf(k ± - fc 2± ) = T d U^(k 2± ) . (76) 
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we can rewrite the first and second terms altogether as 
{p,X,a;q,a,b\Q in ) 



9 e " lPJe " W 7 (2tt) 2 fe?, J (2tt) 2 



x 



(2tt) 2 fc 2 ± 7 (27T) 

^(p, g; fei±, fe 2 ±) + ^(p, g; p ± - fc ± , fc 1± + fc±- pj 



0.6 



(77) 



This is a useful form because, as we will face soon below, the sum of two terms 
in the angle brackets has partial cancellation in some cases. 



5 Rapidity Correlation 



Now that we have worked out the amplitude, its squared quantity immediately 
leads us to an expression for the correlation function shown in Fig. 1 (C) (i.e. 
connected diagram), that is, 



g 2 f d 2 p ± dV d 2 k l± d 2 k' 1± ( P a p (k 1± )pi*(k' 1± )) 



i6tt 2 J (2vr) 2 (2tt) 2 (2tt) 2 (2vr) 2 



k 2 k' z 
K i± K 'i± 



x tr[T ( » -7* ](c/ 6o (fc 2± )C/^'(^ ± )tr[TM 



/ 7^ tr [ T (» •r ( t <) ](c/ 6a (fe 2± )tr[T|C/(fc' ± - k' 2± )T^(k' ± ) 



+ 
+ 



I 



d 2 k ± d 2 k' 



(27r) 2 (27r) 2 "™ 
x (tr[?7(fe ± )71C/t(fe ± - fe 2± ) C/(fc' ± - k' 2± )T^(k' ± ) 



(78) 



where we have suppressed the momentum arguments for = T^(p, q; ki±, fc 2 j_), 



T tM^ = T-W a- k' 



■or 



(>) 



T, 



(<) 



(>) 



"1_L> ^jj) J (<) 



= T^(p,q;p ± - k ± , k 1± + k±- p ± ), and 



T}<\(p,q',P± — k' ± , k[ ± + pj_). We note that k 2 ± and k' 2± are 



defined through k 2 ±= k' 1± + k 2± — p± + q±- These are our final results. 

We present the results in a somewhat more intuitive form in terms of the 
gluon distribution functions in Appendix A. The evaluation for the Wilson 
line average is supplemented in Appendix B. 

In the rest of this section we further explore the rapidity correlation out 
of Eq. (78) in some special situations. First, let us consider the remaining 
correlation which is rapidity independent in the long range. We assumed 
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I Ay I < l/a s to avoid the quantum evolution effect, as we mentioned in Intro- 
duction. Now we shall consider the correlation in the rapidity region satisfying 
1 <C | Ay | < l/a s , so that the eikonal approximation can work for emitted glu- 
ons. Then, we set y 2 1 and y 1 <C — 1 and regard q + and p~ as much greater 
than other scales. Then, q + is so large that the gluon with q + cannot prop- 
agate longitudinally due to suppression by l/q + in Eq. (5). As a result, the 
longitudinal component of the gauge field A" is to be neglected in the eikonal 
limit. 

In this approximation it is immediate to see; 
(p + q) 2 « 2q + p~ , 

C-(p + q,k 1± ,k 2± )&0, C\p + q,k 1± ,k 2± )^-2k\ (79) 
r^(p, q) ~ q + (-9+\9^ + Ig+nQvx) , 

where we have dropped two terms from T^ip, q) because no finite contribu- 
tion arises whenever either of /x, A, or v takes +. After some calculations with 
using the on-shell conditions we can find and to be rewritten as 



^ ) (p, ? ;fe 1± ,fe 2± ) « - V = ^f( p± )f(k 1± ), (80) 



A(p-k) i k{ Aixk 



T l ^(p, q; p ± - k ± , k 1± + k ± - P± ) » k [y = ^^f(P±- *l) P(k 1± ) , 

P± ^ ^ (81) 

where we define /' following the notation in Ref. [10] by 

V n p± 

We note that / i (fcij_)p p (fcij_) is nothing but the perturbative Weizsacker- 
Williams gluon field. From a different point of view, the above together with 
e\ X \p) can be interpreted as coming from the amplitude of the soft gluon 
emission in the color dipole picture. 



After all, the total amplitude becomes as simple as 

(P,ha;q,j,b\a m ) = -— [ f J (k 1± ) p c p (k 1± ) 

g J (Iny 

d 2 k 



(2nf 



f\p ± )T d A U dc {k 2i _) 

(83) 

f(p ± - k ± )U(k ± )T c A U^(k ± ~ k 2± ) 



This production amplitude given in coordinate space exactly coincides with 
the third term of Eq. (5.8) in Ref. [10] except for the normalization. We should 
remark that Eq. (83) is independent of the rapidity difference Ay as we neglect 
e~ Ay under the condition that only q + and p~ are dominating. 
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Next, we will consider a different situation with the rapidity dependence han- 
dled easily. Let p ± and q ± be much bigger than ki±, k 2 ±, and k±. This means 
that p ± and q ± belong to the hard scale as compared to the saturation scale 
which characterizes the gluon distribution inside the target nucleus. We note 
that this situation is the case if we are interested in the di-jet correlations. We 
shall further limit our calculation here to the back-to-back kinematics, i.e. 

p ± + q ± = k 1± +k 2± = ± , (84) 

for simplicity. 

As a matter of fact, the back-to-back kinematics significantly reduces com- 
plexity enabling us to perform explicit calculations with a reasonable effort. 
Without any expansion we can express T(>) and T(<) as 

T^ ) {p,q;k 1± ,k 2 ±) 

~k\ ± tanh(Ay/2) + 2k 1± - p ± }g ij + 2k{p\l + e~ A f ) + 2*V(1 + e Ay ) 



p 2 j_{l + cosh Ay) 

(85) 

^ ( <)(P,?;P±- k ± ,k 1± + k ± -p ± ) 

fen.- (2p ± - 2k ± - k 1± )g ij + 2k{ (p l - k l ) (1+e" Ay ) + (p 7 - - fc? ) (1+ e Ay ) 
~ (p ± - fc ± ) 2 (l + e- A ^)/2 + (p ± - fc_L- fe 1± ) 2 (l + e A ^)/2 ' 

(86) 

Here it is easy to make sure that the above two expressions are to be trans- 
muted to Eqs. (80) and (81) in the limit of Ay — > oo. It is also manifest that 
T(>) + T(<) is vanishing up to the first order in terms of k\± and Let us 
henceforth look more into the rapidity dependence to grasp a feeling from the 
above results. 

We expand Eqs. (85) and (86) in terms of ki± and k± keeping the quadratic 
order. After tedious but straightforward calculations we can reach; 

T ij_ T H _i_ rpij „ 2k 1± -k ±g v + 2fcjF(l + e- A ^) + 2fcjF(l + e Ay ) 

1 + pi(l + coshAy) ' (87j 

where we define 

k± = k±- 2(p ± ■ k L )p ± , k ± = k ± + 1+ 1 e _ Ay ki± , (88) 

and pj_ = Pj_/|pj_|. The system having translational symmetry in the trans- 
verse plane constrains as k' 1± = k±±, but not for k± and k' ± . The squared 
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quantity of the above matrix takes a relatively simple form of 

8 (k 1± - k±)(k 1± - k' ± ) + 2k 2 1± (k±- fc' ± ) cosh Ay (1 + cosh Ay) 



(pi) 2 (l + coshAy) 2 

(89) 

Here the angular average with respect to p ± can take us to make further 
simplification, so that we finally get 



[ J " (Pi) 2 I 1 



2 cosh Ay" 
- cosh Ay 



(90) 



Here we should not forget that k± and k ± have the rapidity dependence too. 
By shifting k± — > k± — k\±/2 and k' ± — > — k\±/2, the final expression for 
the back-to-back two-gluon production from the connected diagram is 



diV„ 



d\p ± \d\q ± \d9dy 1 dy 2 



4a s 



q ± =—p ± l conn 



nRl (Ar c 2 -l)(2vr) 5 pi J (2tt) 



x 



/ 



d 2 fc d 2 k' r 



(2tt) 2 (2tt)^ 



1 2 

fe ± -fc' ± + -(tanh(Ay/2)) fc 



2 



1 + 2 cosh Ay" 
1 + cosh Ay 



(91) 



x (te[U{k 1 _-k 1 j2)TlU\k 1 +k 11 ]2) U(k' i _+k 1 j2)TlU\k' ir k 1 j2)]) , 



where 9 is the angle between q ± and and we have used the gluon distri- 
bution function defined in Eq. (A.l) in Appendix A. The transverse area ttR 2 
originates from our definition of the gluon distribution function. 

The expectation value of the Wilson line correlator includes contributions 
from all twists, which is numerically evaluated. If we make a leading twist 
approximation, we can analytically evaluate the Wilson line correlation, which 
is simply obtained as 



ti[U(k L -k l j2)TlU\k r ±k 1 j2)U(k' ± +k 1 j2)TlU\k' 1 _-k 1 j2)} 
4a s (p A (k 1± ) 



k 2 



1_L 



k 



k 



^) 5 ( 2 )(fc' ± _ 
■-)) 

2' v "^' 2 /_ V "- L 2 ' 2' V "- L ' 2 /_ 2 >> ' 



h^k ±+ ^)S^ ± -^) - >(fc' ±+ ^)^( fcj 



k 



k±± , 
2 ' 

(92) 



where we have defined y?A to mean the gluon distribution function of the 
nucleus. In this approximation, we obtain the correlation function as 



di\L 



i d|p ± |d|q ± |d6'dyidy2 



1 



4a 2 iV 2 



q ± =—p ± l conn 



nRl(N^-l)(2^plJ (2k) 



I 



d 2 k 



i_l 



x 



M k i±) ( PA(k 1± )h+ (tanh(Ay/2))' 



1 + 2 cosh Ay" 
1 + cosh Ay 



(93) 
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2.5 




J4 




Fig. 5. Dependence on Ay in the leading-twist approximation. 

Although the expression in Eq. (93) has a /c_i_-factorized form, we cannot write 
the two-gluon production in general into such a /c^-factorized expression even 
in the pA case [9,12]. In this back-to-back kinematics, the rapidity dependence 
is simply factorized as 



dN B 



d\p ± \d\q ± \d9dy 1 dy 2 
4 



-px ' conn 



X 
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3 + (tanh(Ay/2)) 



d\p ± \d\q ± \d6dy 1 dy 2 

1 + 2 cosh Ay V 
1 + cosh Ay J 



-p ± ,Ay=0/ conn 

(94) 



The correlation function is a monotonously increasing function of Ay. In Fig. 5, 
we show the rapidity correlation associated with the two-gluon production 
in the back-to-back case. It is remarkable that the correlation is enhanced 
as the rapidity difference increases, which is valid as long as Ay < l/a s . 
At asymptotically large rapidity difference with Ay — > oo, if the quantum 
evolution between Ay is neglected, the correlation function exactly reproduces 
the result in the eikonal limit, which is in agreement with Ref. [5]. 



6 Conclusions 



We explicitly calculated the two-gluon correlation in the pA collision using 
the framework of the Color Glass Condensate and the LSZ reduction formula. 
We specifically chose a connected diagram, which is the leading order of the 
coupling g and the weak color source, and could be the only contribution to 
produce the intrinsic dependence on the longitudinal or rapidity separation 
between emitted gluons. To accomplish the calculation we wrote down the 
classical gluon field associated with the dilute and dense color sources and put 
it to the background propagator. 
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Our final results encompass the rapidity dependence that is not taken into 
account in the eikonal limit, but they are too complicated to perceive the 
qualitative behavior as a function of the rapidity separation. We have to wait 
for numerical evaluation, which we plan to do as a natural extension of this 
work, in order to extract phenomeno logical physics implication from our ana- 
lytical expressions. 

In this paper we proceeded to analyzing the rapidity dependence in the special 
cases. First, we confirmed that our expressions are certainly reduced to the 
known expression given in Ref. [10] in the eikonal limit. Next, we relaxed 
this limit so as to allow for the rapidity separation under the condition that 
the momenta carried by emitted gluons are hard as compared to the gluon 
distribution in the projectile and target hadrons. In this case we found a fairly 
concise formula for the two-gluon production in the back-to-back kinematics. 
Our results imply that the dependence on the rapidity separation is only 
moderate and there arise no exponentially damping correlations. 

In fact, our finding is the other way around; the connected diagram enhances 
the long-range correlations in rapidity. That is, the contribution to the multi- 
plicity from the connected diagram gets larger with increasing rapidity separa- 
tion. Although this might sound curious at first, this could be possible because 
we consider two gluons split from one gluon. The long-range correlations in 
this case simply signifies that more gluons come out with larger longitudinal 
angle. In fact, the three-point vertex in QCD contains p + /q + = e Ay which 
surpasses the energy denominator. We remind that our tree-level description 
breaks down eventually for Ay > l/a s due to quantum effects. 

It is quite interesting to see such long-range correlations in rapidity and this 
observation seems to be consistent with the ridge phenomenon observed at 
RHIC [27,28,29]. We cannot, of course, make any convincing statement about 
the near-side jet correlations from our rough analyses limited to the back- 
to-back kinetics. It is indispensable to perform the numerical integration to 
quantify the two-gluon production. In the near future we are planning to 
include the disconnected diagram, which is much simpler than the connected 
one we manipulated here, to discuss not only the rapidity dependence but also 
the azimuthal angle dependence and their interplay in this framework [42]. 
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A Appendix - Gluon Distributions 



To make the final expression (78) more comprehensible in an intuitive way 
let us make use of the notation by the parton distribution functions along 
the similar line as Ref. [9]. We first define the proton unintegrated gluon 
distribution; 



5 aa ' k\ L 



d 2 X 



3 i(fcix-fci ± )-Xj 



dy9 p (fc 1± |x ± ) (A.l) 



d 2 X 



where we have assumed k\± ~ k' 1± . The above equality could be understood as 
the definition for <£ p (k±\X ±). It is a straightforward generalization to proceed 
to the definition for the two-point correlator, that is, 



JV c (tr 



5 aa \U ha {k 2L )U Wb \k' 2i _)ti[T b p Tl 

' " d 2 X, 



U{k 2i _)U\k[ 



2L> 



7lk 2 



2± 



(A.2) 



We should note that b is the impact parameter that specifies the relative 
separation from the proton gluon distribution in the transverse plane. The 
three-point correlator goes on to be 



itk 2 21 _ 



U(k ± )T^(k ± -k 2± )T^]u^ b \k' 2± )) 

gg ' g (k ± ,k 2± -k ± ;k 21 \X 



I 



d 2 X 



j_e 



±- 



b) (A.3) 



d 2 X 



The four-point correlator is 



5 aa '(tr[u(k ± )T^(k ± - k 2± )U(k' ± - k' 2± )T£u\k' ± ) 

g 2 N- 



7Tk 2 2± 



1 



d 2 X i e i ( fc 2±-fc' 2± )-x J 



d v T gg (k ± ,k 2± -k r ,k' ± ,k 2± - 
d 2 X, 



k' ± \X ± -b) 
(A.4) 



Here we assume that the distribution functions are approximately homoge- 
neous in transverse space, which is often the case in the McLerran-Venugopalan 
model treatment based on the Gaussian distribution for the color source. Then, 
the integration over the transverse plane gives rise to the momentum delta 
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functions, which simplifies the final results to be 

/ dN g \ 1 a 2 s N 2 r d 2 p ± d 2 q ± d 2 k 1± <p p (k 1± ) 

\d yi dy 2 / conn ^(^2-1)^7(2^)2(2^(2^ k 2 1± k 2 2± 

x |tr[T (>) ■ T ( t >} ] <//(fc 2± ) + J |^ tr[T (<) ■ T^] ^(fe ± , fc 2± - fc ± ; fc 2± ) 
+ / |^ tr[T (>) • T ( t <} ] ^(fe ± , fc 2± - fc ± ; k 2± ) (A.5) 

+ / (2^ t^) 2 tr[T(<) ' T( ^ )] ^ g ' gg(fe - L> k2± ~ k±] ^' fc2± ~ fc/±) } ' 

where _R P is the transverse radius of the proton. It is necessary to calculate 
the gluon distribution functions for numerical evaluation and we can utilize 
the McLerran-Venugopalan model. Appendix B is devoted to the explanation 
of how to compute the average of the Wilson lines. 



B Appendix - Average of the Wilson Lines 

In this appendix, we evaluate the Gaussian average for the two-, three- and 
four-point functions in terms of the Wilson line in the adjoint representation. 
They are all necessary to compute the gluon distribution functions in the 
McLerran-Venugopalan model. For this purpose it is enough to calculate the 
four-point function because the two- and three-point functions are readily 
deduced by the appropriate contraction from the four-point function as 

^ c (tr[C/(a!i ± )C/ t (x 2± )]) = (tr[U(x 1A _)nu\x 3± )U(x 3± )TZU\x 2± )}) , (B.l) 

and 

(tr[U(x 1± )nu\x 2± )T b A ]U^(x 3± )) 

= {tr[U(x 1± )T^(x 2± )U(x 3± )T^(x 3± )}) . 

Calculating the generic four-point function in the adjoint representation is a 
tedious task, though it is not impossible. Instead, we thus adopt the large- N c 
expansion here. We can express the adjoint Wilson line using the fundamental 
Wilson lines as 

U(x i± ) Pa = 2tr[TgU F (x t± )T«Ul(x t± )} , (B.3) 

where Up(xi±) is the fundamental Wilson line and Tp is the generator in the 
fundamental representation. Using Eq. (B.3) we can then rewrite the four- 
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point function in the adjoint representation as follows; 
tT[U(x 1± )T£Ui(x 2± )U(x 3± )TZrt(xu.)] 

ArpCL H^CL rpOL\ rpO.2 rpa>% rpOLA 

~ Aaict2 Aa 3 Q 4 F aiai F a^ai F a%az F a^a^ 

x A5wM^- bib T§\ b T§ 3 hb T§\ bi (\{ U Fbiai (x i± )U;- bi -Jx i± 



(B.4) 



In the same way as the present authors did in Ref. [45] we shall introduce the 
"initial state" and the "finial state" as 

(a u a 2 , 03, 04, a ± , a 2 , a 3 , o 4 |initial) = ^T^T^T^^T^^ T«\ 3 - a T«\ 4az 
2(^0104^0201^0302^0403 ^0102 ^0203^0304 ^0401 

^0103^0201 ^0304^0402 ( 5aia2 , 5<22a4 ( 5a3ai ( 5a4a3^ , (B.5) 

(finally, b 2 , b 3 , b 4 , b u b 2 , b 3 , h) = %M feft T F A ^T^ 2ft T F ft ^T F ft M4 



fy>l& 4 "&2&3"&3&2"&401 + J^Oblbl "6262 "6363 "6464 



_ Jj~ (^ b ^ b ^hh^b 4 b 4 + ^6164^6262^6363^461) ' 



(B.6) 



and also define the "Hamiltonian" composed of the free part, H , and the 
interaction V; 



O' 2 

H = ^L(0,0) 



V = 



U=l 

( E N ^ r(^ ± , + T F f T F f r( Xi±: x j± ) 



(B.7) 



where we use the same notation as in Ref. [45] as 



(B.8) 



L(x 1± ,x 2 ±) = g A J d 2 z± G (x 1± -z±) G (x 2 ±-z±) , 
r(x 1± , x 2± ) = 2L(0, 0) - 2L(x 1±: x 2± ) . 



(B.9) 



Here Go(xj_ — z±) is the propagator in two dimensions satisfying dj_Go(x± — 
z±) = 5^ 2 \x± — z±). We note that the saturation scale Q' s has a different 
normalization from Ref. [45], that is, we have defined Q' 2 = 2N 2 /(N 2 — l)Q 2 . 
By means of the bra- and ket-vectors, the initial and final states become 

(initial) = iiV c 2 (|4123) + |2341) - |3142) - |2413>) , 
(final] = iV c 2 (4321| + (1234| - iV c (1324| - iV c (4231| , 
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where we define 

(ai, a 2 , a 3 , a 4 , a u a 2 , a 3 , a 4 \ijkl) = -^5 ai afi a2 a^a 3 aja i a l , (B.10) 

c 

whose norm is normalized properly to unity. Here we note that \ijkiys form a 
complete set of the singlet states out of N C ®N*®N C ®N*®N C ®N*®N C ®N*, 
that is, satisfy H \ijkl) = 0, though \ijkl) is not the orthogonal basis. 

Let us consider the large- N c limit and then one might naively think that the 
leading order of the matrix element would be 0(N*) since both the initial and 
final states are 0(N%). This is, however, not true. In the large- N c limit V is 
diagonal and (final | initial) = 0(N%) [14,45]. Unfortunately, since the leading- 
order contribution vanishes, we need to calculate the next-to-leading order in 
the large N c expansion. 

For the purpose of proceeding to the next-to-leading order, we need to know 
the overlap between different states. These are directly read as 

(4321|4123) = (4321(2341) = — , 

N { (B.ll) 

(4321(3142) = (4321(2413) = — . 

JV C 2 

As we mentioned above, the matrix elements of e~ v are 0(1) for the (vanish- 
ing) diagonal components and are (D(l/N c ) for the off-diagonal components. 
Hence the only first equation in Eq. (B.ll) gives a finite contribution to the 
next-to-leading order, while the second one is the next-to-next-to-leading order 
in this counting. 

As a result of the truncation up to 0(N%) the matrix element is expressed as 
(final|e" y |initial) ~ -^((4321|e _y |4123) + (4321| e -^|2341>) . (B.12) 



Henceforth, let us focus on the evaluation of (4321|e ^(4123) appearing in 
Eq. (B.12). By definition e~ v is expanded as 



-v 



oo / i\n 

= Y,^-V n . (B.13) 



n=0 Ul 



We shall operate V onto (4123) and (4321). Using the diagrammatical tech- 
nique developed in Ref. [45] we can find, 

U|4123) = U n |4123) + ^-U 21 |4321) + • • • , 
U|4321) = U 22 |4321) + --- , 
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where 



v 

Q? 

2 

Q' 2 



r(xi±, x 2 ±) + r(x 2 ±, x 3 ±) + r(x 3± , x 4± ) + r(x 1± , x 4± ) 



V 2 i = -y- [-T(aJij., x A± ) + r(x 2 ±, x 3± ) - r(x 2 ±, x 4± ) - r(x 1± , x 3± ) 



V 22 = ^f- [2r(x 1± , X41) + 2r(x 2± , x 3± ) 



(B.15) 



Here the ellipsis in Eq. (B.14) represents sub-leading states such as 

1 3142), 1 2413), and so on, which are higher order. We can go on this procedure 
to make V n act onto |4123), and then we get the following; 



V n \4123) = |4123) + -rr|4321) V^V^f 1 + V 22 V^~ 2 + ■■■ + V 2 n 2 v ) + 



1 V n — V n 

= 14123)^1 + — 14321)^1-^ 22 



N c < ' V u -V 22 



+ 



(B.16) 



Consequently the necessary matrix element is acquired as 
(4321|e-n4123) = ie-.. +7 L_!iL_ ( 



vu_ e -v„ )+0 (±_y (B.17) 



We take the same path to evaluate (4321 |e y |2341) only to find that 

(4321|e _y |2341) = (4321 |e" K |4123> . (B.18) 
Finally we reach the answer, 



tv[U(x 1± )T^(x 2± )U(x 3± )T^(x 4± ) 



(finial I e~ v (initial) = N c 



e"^ 1 + 



V- 



21 



(e~ Vl1 - e 



Vu -V, 



22 



-V22' 



(B.19) 



Now that we have the four-point function, it is straightforward to evaluate the 
two- and three-point functions by substituting (x 2 ± — > x 3 ±,x 4 ± — > x 2 j_) and 
x 4±— ¥ x 3 ±, respectively. That is, 



N c (tr[U(x 1± )U ji (x2±)]) = N 3 e~ V22 = N 3 e -Q? r (*i±,*2±) ? ( R2 q) 



and 



tr[U{x 1A _)T a A U\x 2A _)Ti]U^ b {x 3A _)) 

= N 3 e~ V22 = N 3 e -iQ' 2 i r ( x i±' x 2±)+r{x 2 ±,x a± )+r{ Xl± , X2± )] _ 



(B.21) 



Although these results for the two- and three-point functions are limited to 
large N c , we can reproduce the correct answer by replacing N 3 by N C (N 2 — 1) 
in view of the results available from Refs. [45,47]. 
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